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For a system of two spatially separated monoatomic graphene layers encapsulated in hexagonal 
boron nitride, we consider the drag effect between charge carriers in the Fermi liquid regime. Com- 
monly, the phenomenon is described in terms of an interlayer Coulomb interaction. We show that 
if an additional electron - electron interaction via exchange of virtual substrate phonons is included 
in the model, the predicted drag resistivity is modified considerably at temperatures above 150 K. 
The anisotropic crystal structure of boron nitride, with strong intralayer and comparatively weak 
interlayer bonds, is found to play an important role in this effect. 



I. INTRODUCTION 

If two systems containing mobile charge carriers are 
spatially separated such that direct charge transfer is 
not possible, but close enough to allow interaction be- 
tween the carriers in different layers, the resulting mo- 
mentum transfer will equalize the drift velocities in both 
systems. This frictional effect was experimentally ob- 
served between (quasi) two-dim ens ional electron gases in 
double quantum well structureJi'^. In most of the theo- 
retical work the interlayer interaction was attributed to 
Coulomb scattering, hence the effect now bears the name 
'Coulomb drag' (see Refs. |3H5]). 

Interest in the subject has been revived recently by 
the experimental progress which made it possible to pre- 
pare two-dimensional electron systems based on mono- 
layer graphene. A considerable number of theoretical 
workgSniSl studied Coulomb drag between massless Dirac 
fermionSjWhich effectively describe the charge carriers in 
graphene^. However, a quantitatively correct explana- 
tion of the experimental data is still lackin^^liSJilor 

In the typical experiment, Coulomb drag is studied by 
driving a constant current I2 through one of the layers 
(the active one, labeled by the index A = 2 in Fig. [ij. 
If no current is allowed to flow in the other (passive, in- 
dex 1) layer, a potential difference Vi builds up there. 
In terms of these two quantities, the drag resistivity 
Pd = {W/L)Vi/l2 serves as a measure of the momen- 
tum transfer between the two layers, where W and L 
are, respectively, the width and the length of the layer. A 
theoretical expression for in second order in the inter- 
layer interactio n can be d erived either using Boltzm ann's 
kinetic equatioiPESHHHl or the Kubo formulcPaiiU. 

In the present work, we focus on the interlayer interac- 
tion responsible for the drag effect in heterostructures 
composed of two graphene monolayers and hexagonal 
boron nitride (hBN), see Fig. [l] The large bandgap in- 
sulator hBN has a layered structure composed of stacked 
hexagonal crystal planes. Recently the material received 
much attention as it allows the construction of graphene 
- hBN devices with, in comparison to the much used 
Si02 substrates, favorable high carrier mobilitied^SHgs] 
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FIG. 1. A sketch of the double layer system under consid- 
eration. The two monoatomic graphene layers (yellow) with 
charge carrier concentration ni, 712 are placed at z = and 
z = d and labelled by the layer index A = 1,2, respectively. 
The surrounding space (regions 1, 11, and III) is filled with 
the insulating material boron nitride with hexagonal struc- 
ture (hBN). 



particular, the Manchester group reported the fabrica- 
tion of devices where a few layer thin hBN crystal, ob- 
tained by exfoliation, is sandwiched between two mono- 
layers of graphene2^H2Sl, if such a structure is used for a 
Coulomb drag experiment, the Dirac fermions in the ac- 
tive and passive layer can exchange momentum not only 
via Coulomb interaction but also by phonon exchange 
through the spacer medium. The effect of a combined 
Coulomb-phonon coupling on the drag resistivity has pre- 
viously only been studied for quasi two dimensional elec- 
tron gases in semiconductor systems^^ 



In the following, we first investigate the effects of the 
anisotropy of hBN, where the bonds in between the 
graphene-like planes are much weaker than the in-plane 
bonds, on the electron - electron interaction via phonon 
exchange. We then show that the inclusion of phonon 
exchange into the description of Coulomb drag can sig- 
nificantly alter the temperature, density and distance de- 
pendence of the predicted value for at temperatures 
above 150 K. 
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II. INTERLAYER INTERACTION 

A. Combined Coulomb — phonon mediated 
interaction 

In a two-layer system as shown in Fig. [l] where the re- 
gions I, II and III are filled with a homogeneous isotropic 
dielectric medium, the Fourier transform of the bare (un- 
screened) Coulomb potential between electrons in layers 
A and A' has the form 



2ev 



(1) 



where q — {qx,qy), Cvac denotes the dielectric constant of 
vacuum and eoo accounts for the high frequency screen- 
ing properties of the medium. Apart from this Coulomb 
interaction, the charge carriers in each graphene layer in- 
teract via a substrate phonon mediated interaction. The 
charge carriers from each layer couple to the long range 
electric fields generated by optically active phonon modes 
in the surrounding material via Frohlich couplinJ^lHIZl 
This remote interaction between carriers in graphene and 
optical phonon modes in a substrate medium was found 
to influence the el ectrical conductivity of graphene on a 
dielectric substrateP12ai22I. 

In Appendix \K\ we show that in an isotropic medium, 
the combined interaction between electrons in layers A 
and A' via the effects of a static Coulomb potential 
and virtual substrate phonon exchange is of the form 
of Eqn. Q, with eoo replaced by the frequency depen- 
dent dielectric function e{uj) of the substrate material 
(see Eqn. ([AT])). 

In the following, we specialize to the anisotropic spacer 
material hBN. From its three acoustic and nine optical 
phonon bands, only those that (via dipole oscillations) 
create long range electric fields couple to the graphene 
electrons^. Given the layered uniaxial crystal struc- 
ture of hBN, these (infrared active) optical modes are 
described by a dielectric tensor of the for nP 



e{uj) = diag [e±{uj), ej_{uj), e\\{uj)] 



(2) 



The resonance frequencies wJ[q and lo^q of the two 
retardecP^ dielectric functions 



e_L,||(w) = ei" + / 



(Wto) 



'^TO , 



(3) 



are the phonon frequencies at the F point for trans- 
verse intraplane shear modes with displacements parallel 
and perpendicular to the c-axis of the crystal (aligned 
with the z direction in Fig. [T]), respectively. We make 
the us ual app roximation of dispersionless optical phonon 
bands.l^^'^^'^ The values for the high frequency dielec- 
tric constants Coo, the oscillator strengths / (related to 
the static, cq, and high frequency dielectric constants, 
/ = Co — £00)1 wto and the damping factors 7 taken from 
Ref.iHare hsted in Table [H 



To obtain the combined Coulomb-phonon interaction 

U^y in the anisotropic medium, we solve Poisson's equa- 
tion 

- V • (e • V(f>) = Pfroc/evac 

with pfreo being the free charge density of a point charge 
— e at the origin. With Eqn. ([2]), Poisson's equation be- 
comes 



dz 



and as uf^ 
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oz 



q,z) \ +q ej_(j}iq,z) 



, d) and U[°^ = U. 
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(4) 



X exp 



-qdil 



A generalization of this result to structures where the 
regions I, II, and III (see Fig. [T]) are filled with different 

insulating materials (or air) is straightforward; u[i^ then 



(0) 



involves different dielectric functions than U. 



B. RPA screened interlayer interaction 

To take into account the screening properties of the 
conduction electrons in the graphene layers themselves, 
we employ the standard procedure of solving the Dyson 
equation for the two-layer system within the random 
phase approximation (RPA) (see Ref. 5). This finally 
yields the dressed interlayer interaction 



^^12(9,'^) 



Ui°,\q,u;) 



(5) 



eRPA('7,t^) ' 

The total screening function for the coupled electron- 
phonon system given by (see Ref. 1311 and Appendix [b]) 

CRPA = (1 - t/i*?Xi)(l - ^i?X2) - U[fu^f X1X2 , (6) 

where X1.2 denotes the (frequency and momentum de- 
pendent) polarizability of the graphene layer^^. 

Fig.[2]shows a density plot of |eRPA((7, w)|, using dimen- 
sionless units x — q/kp and y = Lj/^VFkp), where kp 
is the Fermi momentum. The horizontal dashed green 
lines mark the transverse and longitudinal frequencies 
of the infrared active modes in hBN, connected by the 
Lyddane-Sachs- Teller relation^^ w^q/w^q = eo/eoo- For 
small damping 7 <^ wxoi the real parts of ej^ ||(iuj) are 

close to a pole at uj^q and close to zero at lo^q , re- 
spectively. Near these frequencies, the absolute value 
of the total screening function erpa likewise shows an 
abrupt change from high values (light colors) to almost 
zero (dark colors). In regions where IerpaI is small, the 
red lines Re e^pA = show the coupled plasmon-phonon 
dispersion relation of the two-layer system. 



O.Ob 




1.0 1.5 2.0 

X = q/fcp 



3.0 



FIG. 2. Absolute value of the total screening function erpa 
Eqn. ([g]), with ni = n2 — 0.02 nm~^ and d = 8nm. Vertical 
green lines show the optical resonance frequencies ttixo' '^lO' 
tJxoi £md cjlo of hBN (bottom to top). Red curves mark 
the zeros of Re erpa . The dashed black lines show the line 
y = X and mark the region where Imx = 0. The hybridization 
between phonon and plasmon modes is clear. 
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FIG. 3. Drag resistivity versus temperature for various in- 
terlayer distances, n — 0.02 nm~^. The blue curves show 
|Pd| (Eqn. ([7|) including interaction via phonon exchange 
and Coulomb interaction, the dashed red curves show |pcd| 
(Eqn. ([9|) with Coulomb interaction only, and da shed black 
lines the low-temperature asymptote pco'^ (Eqn. (B2l). The 
lowest pair of curves (d — 8nm) is also plotted on a linear 
scale in Fig. [7| 
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FIG. 4. Drag resistivity versus carrier density for various 
temperatures, d = 8 nm. Colors as in Fig. [3] 
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FIG. 5. Drag resistivity versus layer separation for T = 300 K, 
n — 0.02 nm~^. Colors as in Fig. [s] The low-t emp erature 

"'""■r (dashed black ( " 
verges for large layer separation to Pcrf (dotted green line, 



asymptote PcS (dashed black curve, Eqn. (|B2|), con- 

L line, 

Eqn. (lOl). Note that at this temperature and density pcn'^ 



already differs from the full static calculation, pcD- 
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III. RESULTS FOR THE DRAG RESISTIVITY 

In the following, we assume for the sake of simplic- 
ity the same positive carrier density n (corresponding to 
electron doping) in both layers, such that Ep ^ k^T. In 
particular, we do not address the recently reported drag 
at charge neutrality poinlJ^, which was attributed either 
to contributions from higher order perturbation theory^^ 
or to corre lated density inhomogenities in the graphene 
layerfpl22l. 



FIG. 6. The integral kernel /C of Eqn. ^ with x = q/kp = 1, 
d = Snm, n — 0.02 nm~'^ as a function oi y — uj/vpkp for the 
temperatures 200, 100, 70 K (full curves, from top to bottom). 
The curves have been aligned on the left side by dividing 
with the y — >■ temperature dependence (T/Tp)^. At T = 
100 K (green curve) a peak near the resonance frequency oj\q 
appears, at T = 200 K (blue curve) there is and additional 
second peak near uj^q (see the vertical green lines in Fig. [2J. 
Dashed curves show the integral kernel for pcD, where these 
peaks are absent. 
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The drag resistivity then assumes a negative valu^l^, 
and the first non- vanishing contribution to po obtained 
in perturbation theory is of second order in the dressed 
interlayer interactioipEl. In terms of the variables carrier 
density, layer separation, and temperature, and under 
the assumptions that both layers are with high electron 
doping and T ^ TV*^, it reads (refer to Refs. 1^ and [Till 
for details) 



dx I dy IC{T,d,n) , (7) 

^0 



where ag = e^/(47revaci'_F^) denotes the effective fine 
structure constant in graphene and the integral kernel 

kW..c\Ui2{x,yrx'^^{x,y) 



K. 



sinh (y 



2T ) 



(8) 



The function $, defined in Eqn. (|B1|, is related to the 
nonlinear susceptibility of graphene, and restricts the in- 
tegration range in the a;, y-plane to the region uj < vpq. 

In order to estimate the contribution of phonon ex- 
change to the drag effect, we note that the drag resistiv- 
ity pcD resulting from Coulomb interaction only (which is 
usually taken as a measure for Coulomb drag) is obtained 
by substituting the static value of the electron-electron 
interaction into the integral kernel Eqn. (l8|: 



PCD -PD|y(o)(^ ,,^0) 



(9) 



For low temperatures Ep ksT, the res istiv ity pcD can 
be approximated by Pq^'^ oc of Eqn. (B2) (see Ref. 6 
for a detailed derivation), under the additional condition 
kpd, kpdl e\\ ^ 1 (large layer spacing), this can be further 
approximated to yield^ 



large d 

Pen 



h (4)' C(3) {ksTf 
e2 eJ- 7r28a2 (;j„^)2„3rf4 



(10) 



(Note that in the static limit, one only needs to rescale 
d dy/e±/e^^ and ag — )■ ag/^/eXej[ to take into account 
the anisotropy of hBN.) The full blue curves in Figures 
^]|5]show the absolute value oi pn Eqn. ([t]) for different 
parameters T, n, and d, while \pcd\ is shown by dashed 
red curves, Pq^^ by the dashed black lines in Figures 
and [il and the dotted green line in Fig. [i] shows Ppj^'' . 

AsFigures [3] and |4] show, the contribution of phonon 
mediated interaction to the drag resistivity is vanish- 
ingly small at low temperatures, but becomes notice- 
able for T > 150 K, the effect being more pronounced 
the larger the layer separation. This temperature de- 
pendence is due the factor sinh~^[yTp/(2r)] in the inte- 
gration kernel Eqn. ([8]), which suppresses the integrand 
for values of ?/ > T/Tp. Thus at low temperatures, 
the main contribution to the ^-integration in Eqn. ([t]) 
comes from a frequency range where the dielectric func- 
tions in the integrand are still close to their static val- 
ues. However, the phonon contribution becomes notice- 
able at lower temperatures than one would expect, tak- 
ing into account that the energy of the lowest phonon 
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FIG. 7. Effect of the anisotropy of hBN on the behavior of 
drag with temperature, with n = 0.02 nm~^ and d = 8nm. 
The curves isotropy, ey and isotropy, ti_ were computed as- 
suming that the graphene layers are immersed in an isotropic 
dielectric medium, with dielectric functions given by ey and 
ex, respectively (see Eqn. ([2|). The curve anisotropy was 
computed taking into account the anisotropy of hBN as in 
Eqn. (|4|. Solid curves show po, dashed ones pcD. 



mode hLuj,Q/kB ~ 1100 if. It is also interesting to notice 
that in the range from 100 to 250 K, the drag resistivity 
Pd, including the effect of phonons, is closer to the 
behaviour p^^'^ than the purely Coulomb drag result, 
pcD • The plot of /C as a function of y in Fig. |6] shows 
the origin of the phonon contribution to the integral pjj : 
With rising temperature, peaks near the resonance fre- 
quencies uj\.Q and uij^Q appear in the integrand, which 
enhance the magnitude of p^i. 

Fig.|4]shows that the relative effect of phonon exchange 
on Po is larger for high densities. For high values of n, the 
argument of the dielectric functions e{uj) — e{yvp y^im) 
in Eqn. (|8| reaches the resonance frequency already at 
lower values of y. While pcD decreases rapidly with n due 
to increased screening of the Coulomb interaction, the 
modification of the screening function crpa by phonon 
interaction is seen to counteract this decrease at high 
temperatures. 

Finally, Fig. [7] illustrates the effect of the anisotropy 
in hBN that enters pp) through the electron- electron in- 
teraction Eqn. Q. We compare the drag resistivity in 
hBN with that in an isotropic medium with dielectric 
functions ey and ej_, respectively (see Eqn. The dif- 
ference in magnitude between pp) and pcD is seen to be 
greatest in the anisotropic case, where both in-plane and 
out-of-plane phonon modes contribute to the interlayer 
interaction. 



IV. SUMMARY AND DISCUSSION 

We showed that including the electron - electron inter- 
action via phonon exchange into the theory of Coulomb 
drag significantly changes the magnitude of the predicted 
drag resistivity in graphene-hBN heterostructures. For 
large layer separations, the deviations become noticeable 
at temperatures higher than 150 K. 
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TABLE I. Parameters for the dielectric function of liBN (see 
Eqn. ^) taken from Ref.H^ 
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4.95 
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3.61 meV 


0.995 meV 




170 meV 


97.4 meV 



The experimental data in Ref. 44 exhibits two resonances, a 
strong and a weaker one, for each direction of the polarization 
of incident light. The weaker ones are attributed to 
missorientation of the polycrystalline samples. 



present work will be able to check our predictions. 
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As the lowest phonon resonance frequency in the spacer 
material hBN corresponds to a temperature of approxi- 
mately HOOK, our result at first sight seems to be at 
odds with the notion that phonon effects should be pro- 
portional to the thermal population factor of the relevant 
modes. This is indeed the case for other transport phe- 
nomena, like the substrate limited electron mobility in 
graphene, where real momentum transfer from an elec- 
tronic state (in graphene) to a phono n mo de (in a di- 
electric substrate material) plays a rol^^Ml]^ ^j^g decay 
rate of the electronic state is then overall proportional to 
the thermal population of the phonon mode. Our sce- 
nario however involves the exchange of virtual phonons 
in a process that is of second order in the interlayer 
interactiorP, and no decay processes into real phonon 
states are relevant for po- We note that in Ref. 13 the 
effect of substrate phonons on Coulomb drag was consid- 
ered for the case where a material described by a uniform 
dielectric function fills what is our region III of Fig.[TJ and 
a deviation from the low-temperature behavior oi po 
was predicted for temperatures roughly an order of mag- 
nitude lower than the phononic resonance frequency of 
the substrate material. 

Up to date, there remains considerable discrepancy 
between experimental data on Coulom b dra g between 
graphene layers embedded in Si02 / AI2 O jiSUS] and hBlNpl 
and the existing theoretical work. For hBN, the reported 
drag resistivities in the Fermi liquid regime are roughly 
a factor of three larger than predicted, and the results of 
the present paper do not change this situation. The ex- 
perimentally reported dependence of po for d — 6nm 
and n = 0.018 nm~^ up to temperatures of 240 K'*^ does 
not disagree with our results presented in Fig. [3) Ac- 
tually, it appears that in the temperature range of 100 
to 250 K the inclusion of phonon mediated interaction 
brings the behaviour of drag closer to the low temper- 
ature behaviour than with static Coulomb interac- 
tion only. Nevertheless, an extension of the experimental 
data shown in Ref. [20! up to room temperature would be 
needed to distinguish clearly between pjj and pcD- The 
use of other substrate materials, such as Si02, should 
not qualitatively alter the results of this paper. We think 
that future experiments with devices as considered in the 



Appendix A: Frohlich electron— phonon coupling and 
phonon mediated electron— electron interaction 

Throughout the present work, we assume the dielec- 
tric properties of hBN layers forming heterostructures as 
shown in Fig. [1] to be the same as for bulk hBN. 

The Frohlich Hamiltonian describing the coupling of 
electrons to a bulk polar longitudinal phonon mode (in 
an isotropic homogeneous dielectric material) is given 



= / d\p{r)^ M{Q)e^'^-^ (aq - al^) 



where p(r) denotes the electron density operator, Oq 
(aq) the creation (annihilation) phonon operator with 
momentum Q = {q^, qy,qz) and the matrix clement reads 



M{Q) 



1 

^00 



1 



(Al) 



with the longitudinal optical phonon frequency wlo • The 
phonon mediated interaction between electrons is given 
by 



MQ,uj)^M{Q)MiQrDi^o{Q,uj), 
with the bare phonon propagator 
Dlo((3,w) = 2wLo/(w' 



(A2) 



(A3) 



We employ t he usual ap proximation of dispcrsionless op- 
tical phononPamiM. 

The bare Coulomb interaction is given by 



(A4) 



where €00 takes into account the high frequency screening 
properties of the medium. 

With Eqns. ([A2|)-([A4| and the Lyddane-Sachs- Teller 
relatior^ cj£q /cj^q = eo/coo , we arrive at the combined 
Coulomb and phonon mediated interaction 



f/(g,w) = i/c(Q) + V'(Q,w) 



(A5) 
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with e{uj) the dielectric function of the medium, see 
Eqn.|3] Since the Frohhch coupUng is derived in a phe- 
nomenological approach based on the dielectric proper- 
ties of the material, the combined Coulomb and phonon 
mediated interaction simply reduces to the Coulomb in- 
teraction screened by e(aj), as it should. 

In a two-layer system as shown in Fig. [T] the Frohlich 
coupling coupling between bulk phonons and 2D elec- 
trons of layer A is given bjESl 



Mx{q,q,) ^i. 



2evac('72 + g2) 



(A6) 



where q ~ {qx, %) is a two-dimensional momentum vec- 
tor. In analogy to the above, we now get for the combined 
Coulomb-phonon interaction in a homogeneous isotropic 
medium 



where the dashed and wiggled lines denote the bare 
Coulomb and phonon interaction, respectively, and the 
full curves electron propagators (see Refs. [51and[5T|). 



The function <I>(a;, y) appearing in Eqn. ([8| read; 



$(a;, y) = $+(a;, y) Q{y - x + 2)Q{x - y) 

+ ^'{x,y)Q{l-y-\l-x\) , (Bl) 



where 



t/jf^°(g,c^) = V3,A'(9)+^AA'(g,w) 



(A7) 



Although it is possible to generalize the Frohlich 
electron -phon on coupling for the case of anisotropic 
materialP^'^ and inhomogeneous layered materials^, 
the easiest way to obtain the effective electron-electron 
interaction, taking into account the phonon mediated in- 
teraction, is by solving Poisson's equation for the elec- 
tric potential created by a point charge in the dielectric 
medium taking into account the frequency dependence of 
its dielectric tensor. 



± cosh 



2±x 



2±x 



'2±x 



1 



y J y \ \ y 

For the low-temperature approximation of po , the fac- 
tor sinh~^[j/rp/(2T)] in the integration kernel Eqn. 
which suppresses the integrand for values of y > T/Tp, 
allows one to expand the remaining integrand to the low- 
est order of y. The y integration can then be performed, 
yielding 



Appendix B: Mathematical details 

The dressed interlayer interaction Eqn. ^ is the solu- 
tion of the coupled set of Dyson equations 



= O ■ 

1 



1 2 




A=l 



lowT 
PCD 



2 /-a 



h 2Tralg{kBT ) 



dx i e 



-2dx 



x^{A-x^) 



[{x + 4:aes)^ - 16algexp[-2dxJ'n:ne^/el)] 



(B2) 



where acs = otg/y (see Ref. '6'for details). 
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